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We discuss coherent transport of Cooper pairs through a Cooper pair shuttle. We analyze both 
the DC and AC Josephson effect in the two limiting cases where the charging energy Ec is either 
much larger or much smaller than the Josephson coupling Ej . In the limit Ej <^ Ec we present the 
detailed behavior of the critical current as a function of the damping rates and the dynamical phases. 
The AC effect in this regime is very sensitive to all dynamical scales present in the problem. The 
effect of fluctuations of the external periodic driving is discussed as well. In the opposite regime the 
system can be mapped onto the quantum kicked rotator, a classically chaotic system. We investigate 
the transport properties also in this regime showing that the underlying classical chaotic dynamics 
emerges as an incoherent transfer of Cooper pairs through the shuttle. For an appropriate choice of 
the parameters the Cooper pair shuttle can exhibit the phenomenon of dynamical localization. We 
discuss in details the properties of the localized regime as a function of the phase difference between 
the superconducting electrodes and the decoherence due to gate voltage fluctuations. Finally we 
point how dynamical localization is reflected in the noise properties of the shuttle. 

PACS numbers: 



I. INTRODUCTION 

Soon after the appearance of the microscopic theory 
of superconductivitjsi, Josephson predicted a remarkable 
manifestation of macroscopic quantum coherence^ by 
showing that two superconducting electrodes connected 
by an insulating barrier can sustain a dissipationless cur- 
rent. Since its discovery, the Josephson effect has had a 
tremendous impact both in pure^i^ and applied physics^. 
One of the most recent and exciting developments in the 
research based on the Josephson effect is probably in the 
implementation of superconducting nanocircuits for solid 
state quantum computation^. 

In nanodcvices a new energy scale appears, the charg- 
ing energy, and new interesting effects show up due to 
the interplay between Josephson coupling and the pres- 
ence of charging. The Josephson coupling, leading to 
phase coherence between the two superconducting elec- 
trodes can be understood in terms of the coherent super- 
position of different charge states. Coulomb blockade^ on 
the other side tends to localize the charge thus destroying 
phase coherence. The simplest example of this interplay 
is provided by behavior of the supercurrent through a 
Superconducting Single Electron Transistor (SSET)^. It 
consists in a small superconducting island connected, by 
tunnel junctions, to two superconducting electrodes. 

Additional features emerge if the SSET is coupled to 
mechanical degrees of freedom thus combining the field 



of single electron effects with the intensively studied area 
of Nano-ElectroMechanical systems^. Among the most 
interesting devices in this area there is the electron shut- 
tle (for a review see Ref. . In its essential realization, 
a shuttle system consists of a small conducting grain, 
in Coulomb blockade regime, oscillating periodically be- 
tween two electrodes (source and drain). The essential 
condition to characterize the shuttling mechanism is that 
the grain must be in contact with a single electrode at 
any time. Following the original proposal of a normal 
shuttle Gorelik et al^ introduce the Cooper pair shuttle 
where all the device (electrodes and central island) is in 
the superconducting state. Despite the fact that the cen- 
tral island is never connected to the two superconductors 
simultaneously, the Chalmers group has shown that the 
system is still capable to establish a global phase coher- 
ence and hence support a finite Josephson curreni*2ii£. 
The shuttle does not only carry charge, as in the nor- 
mal metal case, but it also establishes phase coherence 
between the superconductors. Differently from the nor- 
mal metal case, the Cooper pair shuttle does not need 
a moving island, it is just a SSET with time dependent 
Josephson couplings and therefore it can be realized in 
the standard SSET with time-dependent fluxes^i. The 
properties of the Cooper pair shuttle crucially depend on 
the decoherence mechanism which is also responsible in 
driving the system toward a steady state. The presence 
of dissipation modifies the current phase relation, but 
does not (in general) destroy the Josephson current^iii. 
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The effect of gate voltage fluctuations has been analyzed 
in Ref. [ill where it has been shown that decoherence can 
even enhance the Josephson current. Additional work 
on the Cooper pair shuttle considered the full counting 
statistics of Cooper pair shuttlingi^, and the possibility 
of observing quantum chaotic dynamics^^. 

The present paper extends our previous works on the 
subjec t^^i^^ . In addition to our previous results we pro- 
vide details of the derivation of DC Josephson current in 
the limit Ej/Ec <C 1, and consider several extensions 
which are important for a connection with possible ex- 
periments. We also analyze the effect of fluctuations in 
the external driving and the effect of an external voltage. 
Moreover we discuss the AC effect and study the inter- 
play of various times scales on the spectrum of the AC 
current. In the opposite limit Ej/Ec ^ 1, which has not 
been discussed in the literature so far, we present analyt- 
ical and numerical results on the dynamical localization 
and discuss its signatures on the Josephson current fluc- 
tuations. 

The paper is organized as follows. In Section |TT] we 
present the model of the Cooper pair shuttle. In Sec- 
tion lllll we analyze the transport properties of the shuttle 
in the Coulomb blockade regime. We present the details 
of the formalism to determine the steady state density 
matrix, and to derive the steady state Josephson current, 
whose physical features arc discussed in Subsection llll CI 
Subsections IIIIDI and IIIIEI are dedicated respectively to 
the effect of fluctuations of the external driving on the 
Josephson current and to the effect of an applied voltage 
bias, the AC Josephson effect. The chaotic regime of the 
Cooper pair shuttle is the subject of Section HVl We dis- 
cuss the dynamics of the charge in the central island, in 
Subsection IIV Al A new feature which appear as com- 
pared to the kicked rotator is an extra phase shift during 
the kicks which is due to the superconducting phase dif- 
ference of the electrodes. This phase shift plays a key role 
since it is responsible for time reversal symmetry break- 
ing whose consequences for the dynamics are investigated 
in Subsection IIV Bl In Sections IIV CI we develop the nec- 
essary formalism to calculate the full counting statistics 
in the chaotic regime. The concluding remarks are pre- 
sented in Section [V] Several technical details are given 
in the Appendices El El and O Throughout the paper 
kB = 1. 



II. THE MODEL 

The Cooper pair shuttle is schematically shown in 
Fig. [TJ It consists of a central island connected to two 
superconducting electrodes and capacitively coupled to a 
gate voltage. The superconducting island is small enough 
such that charging effects have to be included. The 
two leads are macroscopic and their phases (f>L,ii. can be 
treated as classical variables. The couplings of the island 
to the leads are time-dependent. This time dependence 
is given by external means and can be achieved either by 
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Figure 1: (Color online). Upper left panel. Schematic repre- 
sentation of the system described by the Hamiltonian Eq. ([1} . 
It consists of a Cooper pair box coupled through externally 
switched Josephson junctions to phase biased superconduc- 
tors. Upper right panel. Time dependence of the left and 
right Josephson energies within a single period. Lower panel. 
Sketch of the Cooper pair shuttle's cycle. The three intervals 
L, C and R, within the period T — II + + tn + t^, cor- 
respond to the situations: (L) E/^\t) = El, E/^-^t) = 
(interaction time at left lead); (C) E/^\t) = 0, E/'^^t) = 
(free evolution time in forward and backward directions); (R) 
Ej'-^\t) = 0, E/^'>{t) = Er (interaction time at right lead). 



making the island to move or by tuning in time magnetic 
fluxes and gates. This has to be contrasted with the case 
of single electron shuttle where for the implementation of 
the shuttle a mechanical moving island is necessarj*i^ii^. 
The system is described by the following Hamiltonian 



Ho = Ec{t)[n - ng{t)\' 



E 

b=L,R 



Ef\t)cosiip-cj)b). (1) 



In Eq.(IT|) h is the number of excess Cooper pairs in the 
central island and ip is its conjugate phase, [n,0\ = —i. 
The charging energy is given by Ec{t) = (2e)^/2Cs(0 
with Cs(i) = Cg{t) + CL{t) + CR{t) the total capaci- 
tance of the SSET {C^/ji/g are the various capacitances 

indicated in Fig. [T]), Ej^'^\t) are the Josephson cou- 
plings to the left or right lead respectively, and ng{t) = 
Cg{t)Vg{t)/2e is the gate charge which can be tuned via 
the gate voltage Vg. 

The Hamiltonian in Eq. ([T|) is nothing else than a SSET 
with an external drive contained in the time dependence 
of the Josephson energies and of the gate voltage. If the 
time dependence of the coefficients is neglected the sys- 
tem is a SSET whose physics is known both in the case 
of macroscopic junctions (Ej 3> Ec) and in the pres- 
ence of charging effects {Ej <C Ec)r- By introducing a 
time dependence of the coefficients, it is possible to ex- 
plore different regimes. A case of adiabatic change of 
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Ej ' (t) is that of a Cooper pair sluice which has been 
experimentally and theoretically discussed in the litera- 
turei^iii. The shuttling mechanism we are interested in 
is essentially characterized by the sequence of time lapses 
during which the grain exchanges charges with the leads 
and time intervals during which it is isolated from the 
leads. The island is said to be in contact with one of 
the leads when the corresponding Josephson coupling is 
non-zero (with value El, Eh) (configurations L and R 
in Fig[T|). In the intermediate region (configuration C), 

Ef^ it) = Ef^ (t) = 0. Note that both Josephson cou- 
pling are never on at the same time. As in Ref. [l^ we em- 
ploy a sudden approximation (which requires a switching 
time At ^ ?i/£'^(^)) and suppose Ej^'^'\t) to be step 
functions in each region (sec Fig. [1]). In case of a me- 
chanical realization of the Cooper pair shuttle such an 
approximation is well justified due to the rapid decay of 
the Josephson coupling with the distance between the 
grain and the lead. For later convenience we define the 
functions 

eL{t) = e{t)9{tL-t), (2) 

QRit) = e{t-{tL+t^))9{tL+t^+tR-t), (3) 

in order to write the time dependent Josephson energies 
as 

E':l'\t)=Et,Y,^b{t-nT) , be{L,R} (4) 

nGN 

The total capacitance Cs {t) is weakly dependent on time 
at the contact regions^ and therefore we assume it to 
be constant during the intervals L and R (obviously the 
same hold for Ec{t) = Ec)- In the intermediate region 
(C) it is not necessary to specify the exact time depen- 
dence of Ec (t) , as it will be clear in Section IIIII 

In the rest of the paper we study the transport prop- 
erties of the Cooper pair shuttle. The transfer of charge 
is expressed by the presence of a current at left and right 
contacts. The corresponding current operators are, in 
the Schrodinger picture, 

hit) = 2e^smiip~ipL) eUt), (5) 
n. 

init) = 2e^sin{^-^R) enit), (6) 

corresponding to the coherent exchange of Cooper pairs 
between the grain and the left or right lead respectively. 
Due to the periodical external driving, any interaction 
with the external environment leads to a steady state, 
where every observable is periodic. We will essentially 
ignore transient effects and concentrate on the stationary 
values of physical observables. 

A. Cooper pair shuttle with time-dependent fluxes 

Before analyzing in detail the transport properties, we 
discuss a way to realize a Cooper pair shuttle which does 



not require any mechanically moving part. Here the 
time dependence of the Josephson couplings and Ug is 
obtained by a time dependent magnetic field and gate 
voltage, respectively. The setup consists of a more com- 
plicated superconducting nanocircuit in a uniform mag- 
netic field as sketched in Figl^l By substituting each 
Josephson junction by SQUIDs, it is possible to control 
the E^j \t) by tuning the applied magnetic field piercing 
the loop. The presence of three type of loops with dif- 
ferent area, A^, Ar, Ac allows to achieve independently 
the three cases, where one of the two E/s is zero (regions 
L,R) or both of them are zero (region C), by means of a 
uniform magnetic field. If the applied field is such that 
a half-fiux quantum pierces the areas Al,Aji or Ac, the 
Josephson couplings will be those of regions R,L and C, 
respectively and the Hamiltonian of the system can be ex- 
actly mapped onto that of Eq.([T|). Moreover, by choosing 
the ratios Ac/Ar = 0.146, and Ac/Al = 0.292 the two 
Josephson coupling are equal, E^ = Er = Ej. This im- 
plementation has several advantages. It allows to control 
the coupling with the environment by simply varying the 
time dependence of the applied magnetic field. The time 
scale for the variation of the magnetic field should be con- 
trolled at the same level as it is done in the implementa- 
tion of Josephson nanocircuits for quantum computation 
(see Ref. [H for an extensive discussion). 

For a quantitative comparison with the results de- 
scribed here, the magnetic field should vary on a time 
scale shorter than h/Ej, typically a few nanoseconds 
with the parameters of Ref. llSl . This is possible with 
present day technologjii^. At a qualitative level the fea- 
tures of the Josephson current presented in this paper 
do not rely on the step-change approximation of the 
Josephson couplings. These effects are observable even if 
the magnetic field changes on time-scales comparable or 
slower than Ej. The only strict requirement is that only 
one Josephson coupling at the time is switched on. 



III. COULOMB BLOCKADE REGIME 

We first consider the system when E^, Er <gi Ec, i-e. 
in the Coulomb blockade regime. In addition, the gate 
voltage is chosen so that < ng{t) > 1/2. Namely 
Ugit) = 1/2 when the system is in contact with one of 
the lead and ng{t) = const. G (0, 1/2) in the remaining 
time of the cycle. Our choice (the same of Ref. results 
in having exact charge degeneracy during the Joseph- 
son contacts, then enhancing charge transfer. A differ- 
ent condition, of easier experimental realization, in which 
ng{t) = const, is discussed in Appendix [Al In this limit 
one can restrict the Hilbert space of the system to the 
one spanned by the two charge states {\n = 0) , \n = 1)}. 
The Hamiltonian of the system restricted to the two di- 
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Figure 2: (Color online). The setup for the implementation 
of the shuttle process by means of a time-dependent magnetic 
field. The inset shows the time dependence of the applied field 
(in unity of -Bo = 'l?o/(2Ac),$o is the flux quantum) in order 
to realize Cooper pair shuttling. The different loop areas can 
be chosen in order to obtain El = Er. 

mensional vector space, reads 

b=L.R 

(7) 

where we used the 2x2 Pauli matrices ai {i = x, y, z) 
with the standard notation <t± = [a^ ± ay)/2. 

In order to evaluate the current, Eqs. (|5l6p , or the av- 
erage value of any observable, we need to compute the 
reduced density matrix of the central island p{t). The 
steady state density matrix depends on the specific de- 
coherence mechanism. The main source of decoherence 
in the Cooper pair shuttle is due to gate voltage fluctua- 
tions, either induced by the electromagnetic environment 
or by background charges. 

A. Classical noise 

At a classical level voltage fluctuations can be included 
by adding a classical stochastic term to ng{t). The 
Hamiltonian in Eq. ([7]) is modified by the presence of 
the extra term, 

H^Ho+ m'y^ , (8) 
where ^{t) has a white noise spectrum 

(C(i))stoc - 

where 7 is the inverse decoherence time. If we neglected 
the fluctuations, the time evolution of the system would 



be fully coherent. By including fluctuations, the shuttle 
will be described by 2 x 2 density matrix that obeys the 
following Bloch equation: 

§f =-^(i?oWp-/5^oW) -27(/5-a,/}a,) . (9) 

The only stationary solution of this equation is trivial: 
p = 1/2, this corresponds to the absence of any aver- 
age superconducting current. This is a combined effect 
of the decoherence term and Josephson coupling. In the 
absence of Josephson coupling, voltage fluctuations can 
not cause transitions between the charge states so no re- 
laxation takes place. With Josephson coupling switched 
on, the voltage fluctuations cause transitions between the 
stationary states separated by energy E^j^'^K Classical 
voltage fluctuations result in equal transition rates with 
increasing and decreasing energy. The vanishing of the 
critical current has a simple explanation, the classical 
noise mimics a bath at high temperature^. No coher- 
ence can be established at temperatures much higher 
than the Josephson coupling energy. Nevertheless this 
model is worth to be considered because, as shown in 
Ref. [l3 there is an high temperature regime where the 
average current is zero but still coherence manifests in 
the higher moments of current fluctuations. At low tem- 
peratures, Tfc < Ef"'^'^ the interactions with the bath 
can lead to a density matrix p ^ I and then to a non- 
vanishing supercurrent. 

B. Quantum noise 

In order to analyze the low temperature regime, we 
need to take into account the quantum features of the 
bath. As the most important source of fluctuations in 
the charge regime are gate voltage fluctuation we couple 
the shuttle via the charge operator n to an environment 
described by the Caldeira-Leggett model^S, 

Hint = nO + Hbath = E + "I ) + Hbath ■ (10) 

i 

In Eq. pU]) . Hbath is the bath Hamiltonian, with boson 
annihilation/creation operators of the i— th mode a^, a^, 

Hbath = J2i '^iiO'lo-i + 1/2)- Duc to tlic pcriodicity of the 
external driving the time evolution of the system at long 
time, t ^ T, can be determined by iterating the evo- 
lution of the density matrix p{t) over one single period. 
This evolution can be computed through a linear map 
Mt^t+T defined by 

pit + T)^Mt^t+T[pit)]. (11) 

With the following choice of parametrizing pit) = 
1/2 [I + a ■ r{t)], where i = x,y,z and ri{t) = {(Ji), the 
map in Eq. (|11[) assumes the form of of a general affine 
map for the vector r{t): 

r(t + T) = Mtr{t) + v* , r e 6i(0) C (12) 
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where Bi is the Ball of unitary radius in R^. 
Mt fulfills the property 

|Mtv| < |v| Vve6i(0), 



The matrix approximation. The solution of the master equation in 
the contact region can be obtained in the form 

(13) r{tL) = exp (Gl^l) r(0) - 2^lGI^ ■ [I - exp (Gl^l)] . 



as we will see from its explicit form determined below. 

In the long time limit, the system reaches a periodic 
steady state, 



.it) 



(14) 



if, and only if, det(I — Mt) ^ 0. When this condition is 
not satisfied the external bath introduced in Eq. (fTO|) is 
not effective and the system never loses memory of the 
initial conditions. 

The stationary limit is the fixed point of Ait^t+T^- 
The expression of roo{t), and therefore of Poo{t) = 
1/2 [I + (J • roo(t)], uniquely determines the steady state 
of the system. 

The periodic time dependence of any physical observ- 
able A is given by 



{A{t))^Tr{p^{t)A}, 



(15) 



where the operator A is in the Schrodinger representa- 
tion. 

The assumption of a stepwise varying Hamiltonian 
considerably simplifies the form of the map M.t^t+Ti 
which now can be expressed as a composition of the time 
evolutions of p in the intervals L, C, R (see Fig[T]). In 
each time interval it is straightforward to solve the cor- 
responding master equation for the reduced density ma- 
trijs^. In the portion of the cycle corresponding to the 
island being in contact with the left electrode the master 
equation reads 



r(0 = GL(t)r(t) + 27iWt 
with wj^ = tanh(i?L/7h) (cos0l, sin^L, O) and 



Gz 




(16) 



(17) 



Here, 7^ is the dephasing rate (for this portion of the 
cycle), depending on the temperature of the bath, which 
is taken in thermal equilibrium at temperature T;,. The 
master equation when the island is in contact with the 
right superconducting lead goes along the same lines with 
the substitution L —f R (thus introducing a dephasing 
rate 7/?). Both dephasing rate can be obtained in the 
Born-Markov approximation^, which requires that the 
bath autocorrelation time is the smallest time scale in the 
problem. This treatment is valid provided that Jl{r) 
Th/fi,, i?2,(_R)/S, and that the time interval ^^(fl) is much 
longer than both hTf^^ and hE~^j^y As an example, for 
an Ohmic bath with coupling to the environment a ^ 1, 
one has -il(r) = iTr/2)aE^ii-^coth{E^R-^/2Tb)/m. Gl 
is time independent as a consequence of the Born-Markov 



(18) 

The parameters of the Hamiltonian enter the final results 
only through the combinations = -E'L(/?)^L(i?)/^ 

and ^L(R)tL(R)- Due to the condition ^^n) < E^n^/h 
the parameter ?i7L(fl')/£'L(_R) does not enter the results 
at lowest order. 

During that part of the cycle when the island is dis- 
connected from both electrodes, the situation is simpler. 
Since n is conserved, the evolution can be determined 
exactly 

r(t_ + II) = exp(G^t^) • R{x-,)r{tL) ■ (19) 
In the previous equation we defined 



G_ = 



-7- 






(20) 



and 



^cos(x^) -sin(x^) 0^ 
R{X^) = I sin(x-.) cos(x^) 
1, 



(21) 



where = Jt^^^^ Ec{t){l - 2ng)/h. The rate 7^ 
depends only on the properties of the bath. Its explicit 
time-dependence varies when the time scale is compared 
with the inverse ultraviolet bath mode cut-off, 1/ujc, and 
the inverse bath temperature, h/T)^^^. An expression of 
7^ in terms of bath parameters can be obtained within 
the same Born-Markov approximation discussed above 
in the case of a weakly coupling between the bath and 
the system. It gives 7_f = 2TTaTb/h in which case 7^ is 
independent on time and the decay is purely exponential. 
The same equation holds in the backward free evolution 
time 



r(T) = exp(G^i^) • Rix^)r{T - t^) , 



(22) 



where G^ is defined as G^ in Eq. with the replace- 
ment 7^ 7^. In addition to the dynamical phases 
X-t(^} and 1^(^11^, also the phase difference (j) = 4'l — i'R 
enters in determining the physical observables. The effect 
of damping is characterized by the dimensionless quanti- 
ties lL(R)tL(R), and 7^(^)<_(^). 

From Eqs. pGHlOp it is easy to check that Mt ful- 
fills the property in Eq. [13] except for the following val- 
ues (7L,7fl,7^,7^) = (0,0,0,0) or {jL,7R,dL,dR) = 
(0, 0, kTT/2, /i7r/2)), fc, h integers, when det(I - Mt) ^ 0. 
In these cases, the system keeps memory of its initial con- 
ditions and it never approaches the steady state. This is 
however an artificial situation, because other sources of 
dissipation are present and will drive the system to a 
steady state. 
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Let us comment on the assumption of Heavisidc func- 
tions for Ei^(^ji'j{t) we used to determine the steady state 
density matrix. It defines the simplest model to catch the 
features of the shuttling mechanism, i.e. the existence of 
different regimes during a single period time evolution. 
In fact the precise shape of the Josephson energy pulses 
is not relevant, changing it will change the definition of 
the dynamics phases Xi -eqs. (18), (21)- which en- 
ter as parameters in the results for the density matrix. 
What is indeed neglected in our model is the effect of 
exciting higher energy modes and the effects of gate volt- 
age fluctuations during the switching time. These are 
good approximations for switching times At < fi/Ej and 
At'-fi ^ 1 for any dephasing rate i G {^,^, L, R}. 

C. DC Josephson current 

The asymmetry between emission and absorption of 
quanta from the bath leads to a nontrivial fixed point 
(Eq. 114])) for the map M, thus leading to a non vanishing 
Josephson current through the Cooper pair shuttle. The 
current depends on the quantum dynamical evolution of 
the charge on the island and on the interplay between the 
decohcrence and the periodic driving. If, for example, 

I 



and Mx/2, are obtained from Mq, Vq by the ex- 

change of right and left Josephson contacts and of for- 
ward and backward free evolution time. This means that 
= VMq, Vx/2 = with P acting on the param- 
eters 6li(fl,),X-»(^),7L(i?,)iL(i?),7^(^)^^(^)>L(i?.) as: 

T': (L,^,i?,^)^(i?,^,L,^). (26) 

The expression of the current, which depends on all the 
previous parameters can be obtained analytically from 
Ea. (j23p by explicitly writing Mq and vq in terms of the 
various parameters. The current depends only on the 
phase difference between the two superconductors (j)R — 
(j)L- It is an odd fimction with respect to the action of V 
defined by Eq. From this observation follows that, 

even for = 0, there can be a supercurrcnt between 
the leads as long as the evolution over a cycle is not 
■p-invariant. In this sense the system behaves like a non- 
adiabatic Cooper pair pump. 

The main features of the Josephson current in the 
Cooper pair shuttle have been discussed in Ref. [ll| and 
we recall them here for completeness however providing 
a number of new results and additional details. In the 



the period T is much larger than the inverse dephasing 
rates, the shuttle mechanism is expected to be inefficient 
and the critical current is strongly suppressed. In the 
following we will describe a quite rich scenario, depending 
on the relative value of the various time scales and phase 
shifts. 

As charge is conserved by the coupling to the envi- 
ronment, [n, Hint] = 0, current can flow only trough 
the electrodes. Therefore, in the Heisenberg picture, 
/^(t) -|- /i?(i) + h = 0. By integrating over a period the 
average current reads {Ir = —II = I) 

I = Tr{n(p(ti) - p(0))} = Tr{n(p(r/2) - p(0))} . 

We set the initial time within a period at the beginning 
of contact with the left lead. 

Using the steady state density matrix Eq. (fT4|) wc can 
derive a formal expression for the DC Josephson current 
in the system: 

/ = |,z • [(I - A/o)-Vo - z • (I - 7\/t/2)-Vt/2] (23) 

where z is the unitary vector (0,0, 1)"^, • stands for the 
usual scalar product in M^, and Mt/2: Mq, Vq are 

defined in Eq. ([12]). Their explicit form is 



(24) 
(25) 

I 

case of Ol = Or = 6, = = X, 7l = 1r = U, 
7^ = 7<- = 7c, tL ^ tR ^ tj, = = tc; Fig. [3] 
shows a typical plot of / as a function of 6 and cj). De- 
pending on the value of 9 (a similar behavior is observed 
as a function of x), the critical current can be negative, 
i.e. the system can behave as a 7r-junction. The current 
dependence on the various phases is the result of the in- 
terference between different path corresponding to differ- 
ent time evolutions for the charge states in the grain. By 
changing 7jtj and 7c^Cj certain interference paths are 
suppressed, resulting in a shift of the interference pattern 
and ultimately in a change of the sign of the current, as 
shown in Fig|21 

Another interesting aspect of the Josephson current is 
that it is a non- monotonous function of jjtj, i.e. by 
increasing the damping, the Josephson current can in- 
crease. The behavior as a function of the dephasing 
rates is presented in Fig|31 The presence of a maximum 
Josephson current at a finite value of "fjtj can be un- 
derstood by noting that the current is vanishing in the 
strong and weak damping limits. In both the limits sim- 
ple analytic expressions arc available. 



Mo = exp{G^t^)-R{x^)-exp{GRtR)-exp{G^t^)-R{X-,)-eMGLtL) 
vo = -2exp(G^t^) • R{x^) ■ [irG]^^ ■ (I - exp{GRtR))wR + jl exp{GRtR) ■ exp(G^t^) 

•i?(x^) • Gl' ■ (I - exp(GLii))wi] , 
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Figure 3: (Color online). Left panel. Supercurrent (in units of e/T) as a function of the superconductor phase difference 4> 
and of the phase accumulated during the contact to one of the electrodes 9. The other parameters are fixed as: x ~ Stt/B, 
g-7jtj ^ 3^-4^ g-7ctc ^ 4/5 The plot is obtained for Tf, < Right panel. Average current (Tf, < as a function of 
the dephasing rates, with 4> = — 37r/4, = 77r/10, x = 57r/6. As a function of 7jtj, the supercurrent has a not-monotonous 
behavior. Note the change of sign in the current obtained by varying decoherence rates in each time interval separately. 



If the dephasing is strong, / can be expanded in powers of e ''i-w^lir) .^j^^^ g 7_(^)t^(_) g^^^^^ leading order, 
2e 



T 



tanh (^^^ ^-ijLtL+j^t^} sm{20L) sm{(j) - x^) + tanh [^'^ e-("^«*«+'^^*^) sin(26'fl) sui{(j) + x^) 



(27) 

For simplicity we assume that the Josephson energies at left and right contacts are equal, as well as the contact and 
free evolution time. In this case the previous expression is simplified to 



Istrong ~ ^ tanh ( ^ ) e-(7^*-+7c*c) cos(x^ + x^) sin(2^) sin(0 + (x^ - X-)) 



T 



E, 



n 



(28) 



It is worth to notice the presence of a net DC current even 
in the case of (/) = as argued from general argument 
presented before. The role of breaking the 'P-invariance 
is then played by the difference of the dynamical phases 
accumulated in the forward and backward free evolution 
time intervals, x^ ~ X^- If instead we assume a perfect 
T'-invariance, we recover the known expression for the 
DC current of Ref. [ll| 

Istrong - |tanh(^|^)e-(-^'^*'^+--*-) 

cos(2x) sin(26') sin (/« . (29) 

Strong dephasing is reflected in the simple (i.e. cx sin^) 
current-phase relationship and in the exponential sup- 
pression of the current itself. Strong dephasing, in fact, 
suppresses coherent transport over multiple cycles, thus 
giving a corresponding suppression of higher harmonics 
in the current-phase relationship, i.e. a suppression of 
terms oc sin^""*""^ ((/)), mG N. 

For the sake of simplicity, from now on, we present all 
the result in the case of perfect P-invariance of the time 



evolution of the density matrix in a period. This is not 
a serious limitation for the experimental setups. 

In the opposite limit of weak damping defined by 
Ijtj < 7ctc < 1 

2e T f Ej\ ^jtj (cos (/) -I- cos 2x) tan sin (/) 

J-vipnk ^ — tanh — . 

T \Tb J -fctc 1 + cos cos 2x 

(30) 

The current tends to zero if the coupling with the bath is 
negligible during the contact time. In this case the time 
evolution in the intervals L, R is almost unitary, while, 
in the region C, pure dephasing leads to a suppression 
of the off-diagonal terms of the reduced density matrix 
p{t). As a result, in the stationary limit the system is 
described by a complete mixture with equal weights. At 
the point (7jtj,7cic) = (OjO) our model is not defined 
as discussed at the end of Section Hill The limiting value 
of the current in approaching such point depends on the 
relative strength jjtj ^ "fctc between this two parame- 
ters. 

The current tends to zero in both limiting cases of large 
and small 7jij. Therefore one should expect an optimal 
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coupling to the environment where the Josephson cur- 
rent is maximum. A regime where the crossover between 
the strong and weak damping cases can be described in 

I 



In the hmit of vanishing 7jij, Eq. ((3T|) corresponds to the 
situation discussed in Ref. M- Indeed, both expressions 
are independent of the dcphasing rates. The difference 
in the details of the current-phase(s) relationship are due 
to the different environment. 

In all the three cases presented here, Eqs. ([29l [30l [3T|l . 
the change of sign of the current as a function of the 
phase shifts or x is present. 

D. Effect of driving fluctuations 

The expressions for the current, discussed in the pre- 
vious Section, depend on the specific form of the cou- 
pling between the system and the reservoir. The deco- 
herence we considered so far originates from gate voltage 
fluctuations. In addition in the shuttling mechanism an 
unavoidable coupling to an environment producing fluc- 
tuations in the period and shape of the driving is also 
present. We are therefore interested also in considering 
the effect of fluctuations in the time dependence of the 
external parameters on the Josephson current. Having 
assumed a step-like time dependence of the parameters of 
the Hamiltonian, noise in the external driving consists in 

I 



In the previous equation the expressions Mt{i) and Vf(z) 
are those defined in Eq. (fTTj) . The index in parenthesis 
indicates the explicit dependence of both Mt and vt on 
various tf,(i). We refer to T{k) = J2b^b{k) as the pe- 
riod although the time evolution is no longer periodic 
(before averaging); T{k) is in fact the time the shuttle 
takes to complete a cycle, k labels the number of cycles. 
The first term in the right-hand-side of Eq. vanishes 
in the long-time limit ft, — > oo. Averaging the previous 
expression over the fluctuations of the switching times 
according to Eq. ([551) straightforward, leading to 

(r^(s)), = (I-(A/,)J-^(v,),. (35) 



simple terms is the limit 7c ^ 0, for a fixed value of 0. 
For example, at 6 = t:/A the current reads 



(31) 

I 

fluctuations of the switching times. This means that the 
contact times t^, and the free evolution times t^, 
take, at any cycle i, a random value tfc(i), b = L, <— . 

It is reasonable to assume that the fluctuations of any 
switching time are independent on the others. In terms of 
time intervals tb{i), it follows that the fluctuations Ati,{i) 
around the average value tb, are uncorrelated at different 
periods. Within the same period the fluctuations of any 
two distinct time intervals are also independent. Hence 

{Ah{t)), = (32) 

and 

{{Ata{t)nAh{j)r)t = {{AtaW'),{{Ah{j)r)t m 

a i j y a b. The integer valued arguments of Aifc(-) 
labels the periods and the subscript index runs over the 
set L, R,^,^. The average on the stochastic process 
is defined by ( )j. We will discuss later the distribution 
function of Ati,{i). 

By using the same notation in Eq. pip , we can write 
the evolution of the density matrix after a finite number 
of cycles, h ^ 1, as 



(34) 



Note that if we had considered only the external driv- 
ing fluctuations (i.e. neglecting the effect of gate voltage 
fluctuations) wc would have found, in the steady state, 
Poo{t) (X I. In this case in fact the evolution of r{t) con- 
sists in an alternate sequence of rotation on the Bloch 
sphere around the (1, 0, 0) and (0, 0, 1) axes. Due to un- 
certainty in the rotation angles it is a random walk which 
leads, at long time, to a uniform distribution over the 
Bloch sphere. The nontrivial result in Eq. ((35)) arises 
because of the interplay between the two stochastic pro- 
cesses of gate voltage fluctuations and switching time 
fluctuations. They arc independent because there is a 
time scale separation between these two processes: Cor- 



2e ^^^^ / Ej \ 2e-'^-'*-' [2e-'^^-'*-' cos ^ + (1 + er^'^-'^-' ) cos 2x] sin </> 

Y ) (1 + e-27,7t,7)(l + e-27,7«,; C0S(/)C0s2x + e-47jtj) 



h h-1 h-2h-2 

r{t + Tik)) = n ^M>^Mt) + E n ^■^*(^ + + + vt(ft - 1) . 

k=l A=0 A=0/j=A 

I 
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relations in the quantum bath do occur on a time scale We are interested in averaging the current after the 
Tc ^ T, while the time intervals ti,{i) do not fluctuate system has reached its steady state 
within any single cycle. 



lim 



2e 



Tr 



N 



M 



M 



[pm 



(36) 



In order to proceed further we need to specify the distri- 
bution function P{Att{i)). Let us note that the distribu- 
tion is meaningful only if P{Atb{i)) = for Atb{i) < 0. 
We consider 

P(Aifc(z)) = e{AUii) + T)e{T - Ah{i))/2T , (37) 

with r < tb V6 as a toy model: the underlying physical 
idea is that the switching time can be controlled with 
a precision 2r and the switching can happen with equal 
probability in the interval \ti, ~ T,ti, + t]^. Wc do not 
expect that this simple form of the distribution function 



can determine the quantitative details of the Josephson 
current, rather it can grasp the main features of physi- 
cal effects due to imprecision in controlling the external 
driving. For the sake of concretcness. let us consider the 
limit of strong dephasing (Eq. (|29p ). when the expres- 
sion for the current considerably simplifies. The strong 
dephasing leads to a rapid loss of memory of the initial 
conditions. One can suppose that this occurs after one 
cycle independently on the averaging process. It follows 
that, in Eq. 



2e 



N- 



Tr 



N 



(38) 



Ai depends only on the stochastic parameters of the last (= jth) cycle and Poo(O) on parameters of the j — 1th cycle. 
By considering the expansion of the denominator in the previous equation as 



N 



N 



1/(E Ti^)) = 1/{NT){1 ^tbi^)/NT +...) 



(39) 



4=1 



4=1 b 



together with Eq. ([M)) . Eq. ([5^ reduces to 



(/) - lim — (Tr 

^ N^oo NT ^ 



N 



T(/i)+tt(j) 



-I [p. 



+ 0{1/N) 



(40) 



In the iV — > oo limit, only the first term in the previous equation is non- vanishing. It means that, as a consequence 
of the strong dephasing, the effect of fluctuations in the term l/T in the definition of the current are ineffective. The 
final expression for the average current is 



2e 



-7;tL • ^J^l\ 

e ^ sm — - — ) 
^ /t 



(41) 



which, in the lowest orders in the small parameter rEc/fi, reads 



{I)t 



2e 



1 fEcT 



^strong 



+ 



1 fEcT 



^ X sin(26') sin(2x) - cos(26') cos(2x) 

Ec Ec Ec 



r 



(42) 



with the further condition l/T ^ 7,7 ^ Ej/h <C Ec/ti ^ 7c ^ l/T. The current I strong is defined as 



the current in absence of driving fluctuations and is that 
obtained in Eq. ([^5]) . The leading order correction to 
Istrong docs not modify the functional dependence of the 
current on the dynamical phases: It is a simple renormal- 
ization of the pre-factor, 1 I - [Ecrf / {<oh^) ■ Eq. ((42l) 
shows that higher order corrections can instead lead to a 
modifications of the functional dependence of the current 
on dynamical phases. 

E. AC Josephson effect 

When a bias voltage is applied to a Josephson junc- 
tion it results in an alternating current. This is the AC 
Josephson effect which can be derived from the expres- 
sion for the Josephson current / = "^Ej suy{(I)r — (j)L) 
supplemented by the Josephson relation ^((/>_r — <t>L) = 
^V. In the present problem the alternating current will 
not be simply sinusoidal, it is therefore convenient to 
consider its frequency spectrum defined as 

i{u;) = [ dte-'"^' {I){t). (43) 

(in the simple case of a sinusoidal current it reads /(w) oc 
5{uj-2eV/h) for w > 0). 

The application of a finite voltage bias in the Cooper 
pair shuttle gives rise to a quite rich situation, due to the 
interplay of of the voltage bias effect and the underlying 
periodic motion of the shuttle. The relative magnitude 
of the characteristic frequencies, ~ 2eV/h and ^ l/T, 
determine different regimes that we are going to investi- 
gate. 

We set the electric voltage of the left and right elec- 
trode respectively equal to Vl = —V/2 and Vr = V/2. 
There is no dissipative current through the system as long 
as eV < 2 A. In the hmit Cl, < Cg, and {V/2) < Vg 
the Hamiltonian of the system is still given by Eq. ([1]). 
The effect of the electric potential Vi(fl') in the two leads 
can be included, by means of a gauge transformation, into 
time dependent phases of the condensate wave functions, 
|^L(i?)) ^ e2*^^^(«)*/'* IV'L(fl)), or equivalently, 

4>L(R) 4>L[R) + '2eVL{R)t/h . 

The Hamiltonian describing the effect of voltage bias in 
the Cooper pair shuttle becomes 

Hac = Ec{t)[n-ng{t)f 

- EjeL{t)cos{(p- (pL- eVt/h) 

- EjeR{t)coa{ip~ (f)B^ + eVt/h) . (44) 

We are interested in the frequency dependent current 

i{u;) = 2e^ / dte-"^'QL{t)I{t) 
n- Jm 

= 2e^ e"'*'"^ / ^ ^« e-"^'I{s + kT) (45) 
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Figure 4: (Color online). Absolute value of the Fourier trans- 
form of the Josephson current in units of Jo = 2eEj /h. The 
plot is obtained for tj =tc = T/4, EcT/h = 60, EjT/h = 6, 
7jr = ^cT = 0.001. In the inset we plot the Fourier trans- 
form of the Josephson current in a restricted range of fre- 
quencies for 'yjT — ^cT — 0.001 (black dotted line) and for 
7jT = icT = 0.1 (full red line). 



where T{t) = Tr {sin ((^ — 0^ — eVt/h) p{t)}, according 
to the definition of current in Eq. ([5|). 

We computed numerically the time evolution of the 
density matrix of the island and obtained from it the fre- 
quency spectrum of the Josephson current. As a warm up 
we consider the simplest case which consists in neglecting 
the effect of voltage bias. This is the same case consid- 
ered in previous Sections, in which, however, we analyze 
the instantaneous current rather then the averaged one. 
The results are presented in Fig 2) The spectrum clearly 
signals the periodicity of the time dependence of the cur- 
rent signal and of the presence of the Ql function. In 
the steady state the current is repeats periodically and 
then T{s + kT) = I(s) in Eq. ([l5|) . As a consequence the 
current spectrum presents peaks at integer multiples of 
the frequency w„ = 2t:/T, 

i{Lo) = ^ Ar,S {oj - 27rn/T) (46) 

n 

with 

An ~ I ds exp(— 27rms/T)X(s) . 
Jo 

The form of An is fixed by the expression of the den- 
sity matrix at the fixed point through I(s), I(s) = 
y-Rz{(t)L)T^[s), and with r(s) determined by Eq. (flSl) . The 
signal T[s) consists of damped oscillations at frequency 
Ej/h, T{s) (X exp(— 7js) sin(i?js/fi, + uq). The Fourier 
transform of such a signal determines the characteristic 
features of the current spectrum presented in Fig.|4l Now 
\An\ displays oscillations at frequencies T/tj modulated 
by a power-law decay function (cx l/w, for w ^ Ej/h). 
The nature of the other peaks shown in the inset in 
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Figure 5: (Color online). Absolute value of the Fourier 
transform of the Josephson current (all panels) in units of 
7o = 2eEj/h in the different regimes: eVT/h <^ 1 (upper 
panel), eVT/h ~ 1 (lower left panel), eVT/h > 1 (lower 
right panel). The plots are obtained for tj = tc = T/4, 
EcT/h = 60, EjT/h = 6. The values of V, in units of 
Vo = h/{2eT), and 7j = 7c = 7 for the various plots are 
indicated in the legends. 



Fig. [4] is completely different, in fact they are strongly 
suppressed by the presence of dccohcrcnce. We can de- 
scribe the mechanism that originates these peaks if the 
decoherence is absent. In this case any component of the 
density matrix is an oscillatory function with frequency 
Ej/H, Ec/h at contact or free evolution interval respec- 
tively. The matching between these two different func- 
tions at the boundary between consecutive time intervals 
determines the phase ak in the expression for I(s) after k 
periods. We getl(s) oc sm{Ejs/h+ kT{Ec - Ej)/{2h)). 
Substituting this expression into Eq. (US)) , we find the ex- 
istence of peaks in the frequency spectrum of the current 
at = {2TT/T){k±{Ec - Ej)/4TTh), fc e Z, as presented 
in Fig. m The stronger decoherence is, the more sup- 
pressed such peaks are. 

This picture is modified in presence of a finite voltage 
bias between the two superconductors, as presented in 
Fig. m The features of the frequency spectrum in this 
case depends on whether the condition 2eV ^ h/T or 
2eV ^ h/T is fulfilled. For small V, one can suppose 
that the system reaches a steady state which evolves 
in time only through the time dependence of the pa- 
rameter (j> ^ <f) + 2eVT/h. It means we can replace 
J(s + kT) ~ s\Ti{Ejs/h + 2eVkT/h). In this regime the 
spectrum exhibits a splitting of the frequencies of the 
peaks u = {2'K/T){k ± eVT/{nH)). This is shown in the 
upper panel in Fig. Ofor fc = 1, the peak at lowest fre- 



quency. The effect the peaks at higher frequencies is the 
same, therefore we focus on the first peak at w = 2tt/T 
since it is the most pronounced (the values at the peaks 
An decays with increasing the frequency as in the case 
V ~ 0, see Fig. |4]). Fig. [5] also shows the presence of 
other peaks at integer multiples of 2eV/h which cannot 
be interpreted within the simple model just presented. 
Moreover the width of the peaks is weakly dependent on 
7j(C) ■ The presence of higher harmonics in the spectrum 
becomes evident for h/2eV ~ T (see lower left panel in 
Fig. [5]). The presence of higher harmonics, although we 
do not have a detailed simple explanation, is expected 
once we write the Josephson coupling in the Hamilto- 
nianEq. ^ as Ejcos{(l)-2eV/h)cos{ip+{(j)L+(t)R)/2)). 
This shows that Ej is modulated by an oscillatory time 
dependent factor. Even if further complicated by the 
presence of time order operator in the time evolution 
operator, the current X{s + kT) would present terms cx 
sm{Ej cos(2ey(s -|- Kt)/h)) which generates all the har- 
monics. In the voltage dominating regime 2eV/h ^ 1/T, 
the spectrum does not display new effects apart from the 
described V^-splitting and an enhancement of peaks corre- 
sponding to high harmonics (lower right panel in Fig.[S|). 



IV. THE CHAOTIC REGIME OF COOPER 
PAIR SHUTTLING 

Up to now we have discussed the properties of the 
Cooper pair shuttle in the limit of small Josephson cou- 
pling, Ej <^ Eq- In this Section we consider the oppo- 
site regime, Ej ^ Eq- If the external time-dependent 
driving were absent, such limit would not be of great in- 
terest: It simply corresponds to a SSET in the Josephson 
dominate regime, whose physics is already known^. Due 
to the time dependence of the Josephson couplings in 
the Cooper pair shuttle, instead, there is a time lapse in 
which the Josephson energy is vanishing (sec Fig.[T]), and 
therefore the charging effects still play a leading role. 

We will show that the dynamics of the Cooper pair 
shuttle mimics that of a Quantum Kicked Rotator (QKR) 
with the additional free parameter, <f>. The kicked ro- 
tator is a chaotic system in the classical limit. In the 
quantum regime it presents a variety of interesting phe- 
nomena, including dynamical localization^^. Despite the 
long-standing interest in the QKR, only few proposals 
have been put forward and the only experimental imple- 
mentation so far has been realized with cold atoms ex- 
posed to time-dependent standing waves of ligh t^^i^^'^^ . 
The Cooper pair shuttle can therefore provide a remark- 
able implementation of the QKRiS. by means of super- 
conducting nanocircuits^. This allows us to study the 
effects of dynamical localization on the transport proper- 
ties of a mesoscopic device. So far the effect of dynamical 
localization on the current in a mesoscopic systems has 
been discussed only in Ref. [2^ for a quantum dot under 
AC pumping in which it affects the shape of Coulomb 
blockade peaks. 
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A. From classical to quantum dynamics in the 
chaotic Cooper pair shuttle 

We start our analysis from the Hamiltonian in Eq. ([1]) 
which is vahd irrespective of the relative strength be- 
tween Ec and Ej. For sake of simplicity we assume 
throughout this section , thus rewriting Eq. ([T|) 



"5 

as 



Ho = Ecn^ ~ Ej^[cos{ip- (f)L)e{t-nTy 
+ cos(y. - <j)HMt - {2n + l)r/2)] , 



(47) 



where 8(t) = 6{t)9{tj — t) having assumed the sim- 
ple limit II = tn = tj, = = tc- The dy- 
namics of the Cooper pair shuttle reduces to that of a 
QKR when the contacts times are short enough to ne- 
glect the effects of charging energy. In this case the 
time evolution operator during the Josephson contact is 
exp[— liJjij cos((^ — (f>L.R)/fi\- In the rest of the cycle 
the shuttle evolves according to the charging Hamilto- 
nian only. Let us explain this in some detail. Because we 
are interested in the dynamics of the system at long time 
compared with the period T, any physical observable, as 
already noted in Section lllll is determined by the density 
matrix at integer multiples of the period and therefore by 
the Floquet operator which is the time evolution opera- 
tor over a period. Under the assumption that Ej Ec, 
if Ec cannot induce a significant change of ip during the 
Josephson contact, the Floquet operator becomes 



UiT, 0)^JrVJlV 



(48) 



where J and V are the time evolution operators respec- 
tively during the Josephson contacts and the free evolu- 
tion times 



(49) 



where 



k = Ejtj/h and K = {2Ectc/h){Ejtj/h) . 

The condition that the superconducting phase difference 
is left unchanged at the contacts reads nEctj/H < 1—, 
thus establishing a condition on the maximum number 
of allowed charge states involved in the dynamics. This 
condition is essentially independent on the exact time 
dependence of Ej{t). The Floquet operator in Eq. 
is the same of the Quantum Kicked Rotator (QKR) with 
the additional parameter (jj = (pn — 4>l- This shows that 
the physics of the Cooper pair shuttle, in the limit Ej ^ 
Ec, may reproduces that of the QKR, and, for (f) ^ 0, 
provides an interesting generalization. 

The kicked rotator is the first model in which char- 
acteristic quantum effects of classically chaotic system 
have been observed numerically^^i^. As the parameters 
k, K in Eq. P5|) are varied, the dynamics of the QKR ex- 
hibits several appealing phenomena, including quantum 



ergodicity, quantum resonances and dynamical localiza- 
tioii^^. For a discussion about classically chaotic system 
and various characteristic features of their quantum ver- 
sion we refer to the literature (see Ref. [2^ and references 
therein). Here we note that the exponential localization 
of the wave function due to interference effects is one 
of the most relevant of the mentioned features. The dy- 
namics of the quantum kicked rotator follows the classical 
exponential instability (characterized by a positive Lya- 
punov exponent A) up to the Ehrenfest time Iet^- This 
sets the time scale at which quantum interference effects 
starts to be important leading to weak localization cor- 
rection to the classical behavior—. After a localization 
time t* the classical-like diffusive behavior is suppressed 
by quantum effects^i^. Since typically t* ^ Ie, the 
diffusive behavior is possible also in the absence of expo- 
nential instability. 

The Floquet operator in Eq. can be written 

through the redefinition p — {K/k)n in terms of 

V' = exp{-if/2k) (50) 
■^'l{R) = cxp[-iXcos((^ - 0i(fl))/*] , (51) 

so that U{T, 0) depends only on K, while [p, 0] = 
—iK/k ~ —tk; "fc plays the role of an effective Plank con- 
stant. The classical limit is therefore obtained for k oo 
(or equivalently — > 0) , keeping K = const and the clas- 
sical dynamics depends only on the parameters K and 

In the classical limit the equations of motion of the 
Cooper pair shuttle correspond to a slightly modified 
Chirikov map 



Pt = Pt-i - A'sin{0t_i - [(i + 1) mod 2]0} , . 
Ot = Ot^^+pt, ' ^^^^ 



with 9 = ip — (f>L and where subscripts label the number of 
contacts with the leads (kicks) . The dynamics of a given 
distribution function in the phase space under the action 
of the Chirikov map at = is known: For K > 1, the 
angular variable 6 is uniformly spread over [0, 27r] after 
few kicks; p follows a diffusive behavior. The role of (j) in 
the classical map can be investigated by following Ref. [s^ . 
The idea is to substitute the deterministic description in 
Eq. ([5^ with a probabilistic one where a random term 
66 1 is added to the second equation of ([5^ . We obtain 
a diffusive dynamics for the charge on the central island 
for K > 1 (details of the calculation arc in Appendix |B|) : 

((nt — JT-o)^) 2Dt, with t measured in integers multi- 
ples of the period and where D is the diffusion coefficient 



/r2 

D ^ - 
2 



1 -2cos{2(j))J2{K) 



(53) 



The QKR follows the classical diffusive behavior up to 
the localization time t* . Quantum interference effects, 
as shown in Fig. IIV Al (upper panel), for t > t* , suppress 
this chaotic diffusion: The wave function is exponentially 



, 2 

<An > 




Figure 6; Upper panel. {(An)^) = {{n — (n))'^) as a function 
of time for K = 10, k = 15, and phase difference (from bot- 
tom to top) (^/(27r) = 0,0.05,0.1,0.4, 0.8,0.25. Lower panel. 
Saturation value of {(An)^) oc as a function of (f) for the 
same parameter values as in the upper panel. 



localized in the charge basis, over a localization length £, 
and we have i ^ t* ^ D^. 

The charge fluctuations of the central island freeze in 
time. The localization length can be further tuned by 
changing the phase difference as demonstrated in the 
lower panel of Fig. IIV Al Such quantum effects have not 
yet found a complete analytic explanation. An impor- 
tant step into the problem has recently been achieved by 
Tian, Kamcncv and Larkin^. They calculated the quan- 
tum corrections at time t >tE where such corrections are 
small (but nonanalitic) in "k. From their approach it is 
clear that the presence of 0, by breaking the time rever- 
sal symmetry in the system, docs affect quantitatively 
the quantum corrections to localization. 

B. Time reversal symmetry breaking - COE to 
CUE crossover 

Although the diagrammatic approach discussed in 
Ref. m shades light on the role of the phase difference 
in the system, it cannot give quantitative predictions on 
the fully localized state at long time. We are interested 
in discussing the effect of (j> in the localized state. The 
breaking of TRS by the superconducting phase difference, 
(j), can be seen by direct verification that the Floquct op- 
erator in Eq. (|48p at (/) = is invariant under the action 
of 



t - 


> -t 




n — 


> n , 


(54) 


f - 


> -if 





while such symmetry is broken at 7^ 0. 




Figure 7: (Color online). Left panel. Probability distri- 
bution for the quasi-energies level spacing at = (green 
squares) and at = 7r/2 (red triangles). The other curves are 
obtained from RMT in case of COE (fuU line), "2-folded" 
COE (dotted line) and CUE (dashed line). Right panel. 
\\P4>{s) - -Pcue(s)||^/||-Pcoe2(s) — -Pcue(s)|P as a function 
of 0. II ■ II is the norm we use to characterize the dis- 
tance between curves. Pcoe2(s) and Pcue(s) are the prob- 
ability distribution of level spacing obtained from RMT for 
"2-folded" COE and CUE respectively. is the probability 
distribution function obtained numerically at different values 
of 0. For both panel k = 200 and K = 10. 

In the analogous problem of localization in disordered 
metals, it is known that the TRS breaking results in a 
reduction by a factor 1/2 of the the weak localization cor- 
rections and in a doubling of the localization length (see 
Ref. Issh . The same kind of effects have been observed in 
chaotic systems as well^. In fact, we also observe a dou- 
bling of localization length at (f) ~ Tr/2 with respect to the 
= case as shown in Fig. IIV Al Note that the effect of 
doubling of the localization length has to be added to the 
effect of variation of the localization length through the 
change in the classical diffusion coefficient in Eq. (|53p . 

The effect of TRS breaking can be further investi- 
gated by looking at the distribution of quasi-energies 
spacing. The reason behind this is the conjecture that 
quantum properties of classically chaotic systems are well 
described by the Random Matrix Theory (RMT)^. This 
has been shown to hold for a wide class of chaotic sys- 
tems^ (though exceptions exist^). In the RMT the 
breaking of the time reversal symmetry consists in a 
crossover from the Circular Orthogonal Ensemble (COE) 
to the Circular Unitary Ensemble (CUE) for the Floquet 
operator. 

In order to check this conjecture we look at the level 
spacing probability distribution function P{s) of the 
quasi-energies of the Floquet operator and compare them 
with the predictions of RMT. The results are presented 
in Fig. [7] where our numerical results are shown together 
with the universal (no fit parameters) curves predicted 
by the random matrix theory. 

The plots have been obtained by considering a Floquet 
operator of 2^" levels averaged over 10 realizations cor- 
responding to random values of {(f)^ + 0i?J/2 distributed 
in [0,27r). 

For (j) = the probability distribution function is in 
perfect agreement with that of a folded spectrum^. This 
corresponds to the fact that, at = 0, the Floquet op- 
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erator of our system is U{T,0) = Fi ■ F2, where -Fi, 
F2 are statistically independent Floquet operators of the 
usual QKR. The analysis of the transition between the 
two ensembles is presented in Fig. [3 right panel. We 
have analyzed the crossover between the two ensembles 
by looking at the distance, at a generic (j), between the 
distribution ant the distribution corresponding to the 
Circular Unitary Ensemble, Pcue- The distance between 
the distributions is define by the norm: 



\\f{s)-9{s)\\ 



1/2 



The results are plotted in Fig. [3 right panel. At = 
the distribution coincides with that of a "2-folded" Circu- 
lar Orthogonal Ensemble. At ~ 10"'^ the distribution 
stabilizes to that of a a CUE. The crossover between the 
two ensembles is occurs at ~ 1/%/D = \/2/k ^ 5 10~^. 



C. Transport properties in the chaotic regime 



Up to now we have discussed the chaotic dynamics by 
looking at the fluctuations of the charge on the central is- 
land. The possibility to observe the various fundamental 
aspects of quantum and classical chaotic behavior in the 
Cooper pair shuttle has been discussed in Ref. [H. We 
proceed further in this direction by investigating how the 
chaotic dynamics affects the transport properties of the 
Cooper pair shuttle. The most convenient way to achieve 
our purpose is provided by the determination of the Full 
Counting Statistics (FCS) which was already analyzed in 
the normal shuttle by Pistolesi^ and in the charge regime 
of the Cooper pair shuttle by Romito and Nazarov^iSi. The 
generating function of the current cumulants at left and 
right contacts, defined through 



ar3(A,A,r)|,. 



dti-dtn{l{h)---Iitn) 



is expressed in terms of the Hamiltonian of the Cooper pair shuttle Hq in Eq. ([17|) by 



^XL.Xuir.O) = Texp-'fo ^a,/2,.«/2(.) 4f , i?A,/2,A«/2 = Hoi^L -^H- XlAr -^(f>R- Afl) . 

I 



(55) 

(56) 
(57) 



The necessary steps to calculate the FCS are summarized 
in Appendix [Cl According to Eq. ([57]) . the task is to cal- 
culate the time evolution for p^-^^^-^-^h) which is a modified 
density matrix in which bra and ket evolves according 



to two different Hamiltonians. To this aim we generalize 
the diagrammatic approach presented in Ref. |42| . In the 
basis of charge eigenstates we write 



, n'_ 

G^Zn^,n'_(t) - {n+\mx.,+xM{nUu'X.-,V^n^); (59) 

t is the time counted in integers multiples of the period T and Ux^^\j^ is given in Eq. ((48|) after the replacing 
4>L 4>L — ^L, 4>R ^ 4'R — Xr- In Eq- ([SQ]) we have also added an operator for the free evolution after the last kick 
for computational convenience; it does not affect the physics because it simply means that we look at all observables 
just before a kick, rather than just after it. We now turn to the Wigner representation for the kernel G, 

"+—"-- n^— n'_ 

I 



where n = (n_|_ + n^)/2 and n' ^ (n^ + n'_)/2 are the Wigner representation with the classical equivalent ker- 
"center of mass" coordinates. We compare the kernel in 
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nel for the evolution of a distribution function in the 
phase space. The classical dynamics is unstable in the 
Q direction^^, and therefore we approximate the quan- 
tum kernel G by its average over Q' . As a consequence 
the inverse Fourier transform of Eq. (j60p picks up non- 
vanishing terms only at 71+ = n_ and n'_^_ ~ n'_ . In doing 
so, we have reduced Eq. (j58|) to an evolution equation for 
the populations of the modified density matrix, 

P^^,k^'"\t) - 5]G(^-"-)(n,n';i)pi^--^")(0) , (61) 



with 



{n'\U 



FV), (62) 



and therefore 



3(Ai,Afl,t) - Tr{/^-^«)(i)} 

= Y,G^^^'^^\n,n';t)pn',n'iO), (63) 



We expect that G^^^^-^^^ (n, n'; t) depends only on the 
difference n ~ n' , and we then average out (n + n')/2 as 
follows 

7,(A.,A«), ^ G(^-^-)(n + /,n' + /;i) 



2Ln + l 



(64) 

The key observation to construct the perturbation theory 
is that the average of the kernel can be reabsorbed in an 
average over the free evolution operators by means of the 
relation 



{n + l\V\n' + 1) = {n\W\n') , 



(65) 



(nlWln') 



exp 



-i-k{n - 1)^2] Sn^n' = W{n)6r^^r.' 



(66) 

The kernel G'^^'^'^"^ {n — n'; t) involves products with an 
equal number of W and W\ therefore the average is 
completely defined by 



W{ni)W{n2) . . . W{n,)W{n[)*W{n'2)* . . . W{n'^* 



exp 



(67) 



The potential W of the kicked rotator is not Gaussian dis- 
tributed. However, in the diagrammatic expansion we are 
going to perform, we will restrict only to two point corre- 
lation functions. Therefore the r.h.s. of Eq. (p7)) reduces 
to 6n.n' and the effects of having a non-Gaussian distribu- 
tion don't appear—. Ones the kernel G^^^'^^-'(n, n'; t) in 
Eq. has been replaced by the averaged one which 
depends only on the difference n — n' , the sum over 
n' in the same equation can be performed noting that 
J2n' Pri' ,n' {0) = Tr p{0) = 1. We now turn to phase 
and frequency space by Fourier transforming the kernel 

g' "'\n — n';t) both in n, n' and in t, obtaining 



3(AL,A«,t) = lim / e-*""0^^"''"^(g;c.). (68) 



The kernel 25 is defined through 



f''''^ {9, 9', q, q'; c., c^o) = 6^'"''"^ (q; u;)S{q - q') . 



with 



(69) 



'^,9',q, q';uj,uja) 



\W I 



(+Ai,,+Ah) 



(70) 



9± = 9 ± q/2, e'^ ^ 9' ± q'/2, tj± ± w/2 and 

Uxl.Xr = JrWJlW. The dependence of on 
Ai(fl) is through Jl^h = f,-^kcos{v~('l>L.R+\L,R)) ^ jt^^ = 



gikcos{ip- 



-Xl.r)) 



The diagrammatic perturbation 
theory consists in a series expansion of the r.h.s. of 
Eq. (|70|) in Jj^WJlW . Any term of the expansion cor- 
respond to a diagram which consists of propagators. 



^iu)T —ik cos(e+i-0b — At) ifccos(e— |— 0i,-|-A(,) 



(71) 



and averaged vertices 



w 



-)-))■ 



(72) 



We have absorbed the frequency dependence in the ex- 
pression for the propagator, from which it is evident that 
the final expression for the kernel is independent on ljq. 
The average process is indicated by a dotted line in the 
diagrammatic expression. In any diagrams, correlations 
can take place between any pair of vertices one of which in 
the upper line and the other in the lower line. Any other 
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term involves expressions of the form WW or W^W^ and 
vanishes accordingly to Eq. (|67p . 

For fc 3> 1, diagrams involving crossing correlations 
lines, are parametric small in 1/k, i.e. in "S. This al- 
lows us to approximate, at lowest order in "S, which is 
the classical limit, the kernel with the sum of all ladder 
diagrams, we will refer to as diffuson. 



J 



Classical limit 



The diffuson is defined through the series of ladder 
diagrams 



V 



J" 



w 



W 



(73) 



The computation of the diffuson become simple due to Eq. ([7^ . It allows to factorizc any diagram as products of 



w 



0^+ 'I" Oi_ 



Jo(2/csin(g/2 - Al)) Jo(2/c sin(g/2 - A;^))e 



iuT 



(74) 



then the series is a geometric series and we find 

oo ^ 

^ 1 - e'"^ Jo(2fcsm(q/2 ~ Al)) Jo(2fc sni(g/2 - A_r)) 



i=0 



By replacing ©'^""^"^q, cj) with ^(^^'^"^(g, g'; w) in Eq. (USD, we finally find 

3(AL,Afl,i) = e-6(^-^«)* ^ [Jo(2fcsinAL)Jo(2fcsinAi^)]* 

I 



(75) 



(76) 



It defines all transport properties of the system. In par- 
ticular the probability of charge transfer per cycle at a 
given contact is 



1 f 



SN,2MJll{k) . 



(77) 



Mel 



The properties of Bcssel functions guarantees the nor- 
malization condition "^j^Pn ~ 1- Probabilities are non- 
vanishing only for even numbers of electrons, i.e. Cooper 
pairs. Note also that the very same existence of such 
probabilities means that the transport is completely in- 
coherent and the chaotic dynamics is efficient in destroy- 
ing any features of the coherent superconducting nature 
of the transfer. 

The limit Xl,b. — > of Eq. ((75)) gives us the kernel for 
the evolution of the density matrix in the classical limit. 
In particular, for long time and large n — n' {lu 0, 
<Z^0), 



V{q,q'-u) 



1 



{kqY/2-iujTt 



6{q 



(78) 



which expresses the diffusive behavior in n—n' with diffu- 
sion constant D = fc^/2 in complete agreement with the 
results obtained from the classical map. It correctly re- 
produces the result obtained in Ref . [43. The dependence 
on (j) disappears from the final expression. In fact we 
do not recover the classical effects in Eq. ([5^ because we 
just restrict our considerations to single impurity correla- 
tions, while the dependence on <f> comes from higher time 
correlations (see AppendixlBj). It is possible to reproduce 
the classical corrections to the diffusion constant within 
the diag rammatic approach. This has been performed in 
Ref. for = 0. 



Beyond the classical limit 

The calculation of quantum corrections to the clas- 
sical limit of the FCS in the diagrammatic approach is 
based on writing a Dyson equation for the four point ver- 
tex q''^^''^'^\9,9' ,q,q'-,Lj,ijJo), as discussed by Altland^ 
in the case of Al = Xr = 0. In the limiting case 
Ai = Ajf = Eq. describes the time evolution of 
the density matrix of the system and, due to probability 
conservation, the leading quantum corrections at (/) = 
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are determined only by the most infrared divergent term crossed diagrams and, in the general case Xl, A/; 7^ 0, is 
of the diagrammatic expansion, the cooperon^. defined by the diagrams series 

The cooperon consists of the the sum of maximally 



11 ■ 



+ 



X 



— - 



" e-j - - - . ff 



„ + 



(79) 



The cooperon diagram is topologically equivalent to a diffuson diagram after time reversion of the lower line, 



8_ 



X 



X 



IT f 



f I 1 



w 



f' 



01. 1 3_ 



(80) 






In this form the correlation lines do not cross and it is clear that the contribution of the cooperon series is at same 
order of the diffuson's one in the 1 /"fe expansion. The sum of the diagrams series in Eq. ([50|) leads to 



1 - e*-^ Jo(2fc sin(i5/2 - (Al/2 - Ai?/2))) Jo(2fc sin(i?/2 + (j) - {Xr/2 - Al/2))) 



5(l?-79') 



(81) 



with = e+ + 9'_,^' = e'_^ + e-. 

If Al = Xji = the cooperon in the low frequency, 
small angle limit becomes 

(82) 

At = the cooperon is divergent in the infrared limit, 

^ 0, ^ 0, and the quantum corrections can be com- 
puted. This is not the case a.t (j) ^ (form Eq. ([5^ . 
limtj^o,e^o C(i?, 1?'; w) is finite) where, due to the break- 
ing of time reversal symmetry, the diagrammatic ap- 
proach does not work^. In fact, by Fourier transforming 
the previous equation with respect to time and phase it 
can be easily seen that the diffusive behavior is modi- 
fied by the presence of introducing a relaxation time 
= [(^0)2/2]. At r0 ~ 1 , which is (/> ~ 1//D, the 
cooperon contribution is completely suppressed. 

The described procedure cannot be generalized to in- 
clude non-vanishing counting fields to determine the 
quantum corrections to the classical result for the PCS. 
In fact such a procedure for Xl = Xr = crucially 
relays the possibility of writing a Dyson equation for 
the vertex function and ultimately on the existence of 
a Ward identity corresponding to probability conserva- 
tion^, Tr{p(i)} = 1. In the generalized case, due to pres- 
ence of the counting fields in the kernel G^^'^'^'^\n, n';t), 
2ij.|p(Ai,,AH)| _^ 2^ similar conservation law cannot be 
established when the counting fields are inserted in the 
propagators. Let us note that, even if we had an explicit 
expression for the quantum corrections to the classical 
FCS; it would coincide with the result of Ref. [i^ in the 
limit Xl — Xfi = = 0. In this specific limit the quan- 



tum corrections are shown to be divergent thus imply- 
ing that they cannot be treated as a perturbation of the 
classical limit. Therefore such a perturbative calculation 
cannot give us quantitative information on the quantum 
steady state of the system^, i.e. the dynamical localized 
state. 

Anyway our results for the FCS give us direct predic- 
tions for the current and current noise in the system. 
From Eq. ^ we find {II) = {Ir) = 0. The low fre- 
quency noise at various contact is easily determined 

Sl{oj = 0) = Sr{uj = 0) = fcV4 . (83) 

showing a signature of the diffusion of the wavefunction 
in the charge space. The correlation between left and 
right currents is vanishing: 

/ {lL{0)lR{t))^0. (84) 

JR 

These results are valid in the classical limit "fe ^ 0. Real- 
istic values of the parameter can range from {K ^ 10, k ~ 
15) to {K > 30, k ~ 10^)). In any case with a finite non- 
vanishing i:, after a time t* ~ TK^, the system reaches 
a steady state with a localized wavefunction. In this 
regime we can easily estimate the average value of the 
current and the fluctuations of the current signal. This 
is achieved assuming random uncorrelated phases for the 
localized wave function of the grain: 

^ n— — l 

P{an) = l/(27r) an G [0,27r). (86) 
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It follows that 



hiR) = 0, (87) 

V^=5vlfTT' ^^^^ 

where the " defines (only in Eq. ((87l [88])) the average 
over the random phases of the usual Josephson current 
= 2e|; (sin((^ - (f'L{R))) defined in Eq. dHllSl). Also 
in the localized state the current is suppressed by the 
chaotic dynamics while the typical value of the current 
fluctuations is an indirect measurement of the localiza- 
tion length. 

V. CONCLUSIONS 

The Cooper pair shuttle is a very rich system where 
to study various aspects of coherence quantum dynam- 
ics of driven systems. In this work we concentrated on 
two different regimes where the charging energy is either 
much smaller or much larger than the Josephson coupling 
energy. In both cases we analyzed the transport proper- 
ties as the Josephson current or the current fluctuations. 
In the charge regime we evaluated both the DC and AC 
Josephson effect. In the other regime we analyzed the 
consequences of the underlying classical chaotic dynam- 
ics on the full counting statistics. 

We believe that the possibility to realize periodic driv- 
ing via time-dependent fluxes opens the possibility to 
study the very exciting area of quantum driven systems 
in variety of different situations beyond that considered 
in this paper. 
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Appendix A: EFFECT OF NON-DEGENERACY 
DURING THE CYCLE 

In this Appendix we discuss the Josephson current in 
the Cooper pair shuttle when the gate voltage fluctuates 
around a fixed value away from the degeneracy point even 
during the connection to the electrodes. In Section Hill we 
assumed ng{t) = 1/2 during the Josephson contacts to 
enhance the Cooper pair transfer and ng(t) = const. 7^ 
during the free evolution time. From an experimental 
point of view, however, it would be easier to avoid this 



periodical modification of the gate voltage Vg. It can be 
therefore interesting to have an expression for the DC 
Josephson current in the case of constant gate voltage. 

If ng{t) = const. = 1/2, the Josephson current can 
be obtained from the expressions of section IIIII with the 
replacement = = 0. If ng{t) = const. 7^ 1/2, in- 
stead, the general expression for for the current (Eq. 
still holds, but differences arise in the explicit form of 
matrix Mq and vector vg. During the free evolution 
time intervals the dynamics is unchanged compared to 
the case discussed before. In the L and R regions in- 
stead, we have to include in the Hamiltonian the term 
proportional to the charging energy difference, Ec, be- 
tween the two charge states. T Hamiltonian (system and 
bath), in the basis which diagonalizcs the Cooper pair 
box Hamiltonian now reads 

E 

H = -cr^ + d{cos{2n)a, + sin{2^i)a^) + Hbath , (Al) 
where O is the same of Eq. ^TU\i . and 

E ^ {El + E^f^ , cos{2 fi) ^ Ec/E . 

The Hamiltonian in Eq. (|Aip has been widely studied—. 
In Born-Markov and Rotating Wave approximations the 
time evolution of populations (diagonal terms) in the re- 
duced density matrix and coherences (off-diagonal ones) 
are still independent. The respective decoherence rates 
read 

7P°P- = 27,,sin2(2^), (A2) 
7f--- = 7jsin2(2Ai)+r,;Cos(2M) . (A3) 

We notice that we have to introduce two different de- 
phasing rates, jj and Fj. In our approximation, they 
are 7j = {tt / 2) a E coth{E /Ti,), and Tj = 2TTaTi,, in case 
of weak coupling of the system to the bath (a <C 1). De- 
pending on the relative strength of the two energy scales 
Ec and Ej, we have different effects. If Eq = 0, (cor- 
responding to 2/i = 7r/2) we recover the Hamiltonian of 
the early case described in section IIIII If Ec <C Ej, 
we have corrections of order Ec/Ej in our previous re- 
sults, and we are not interested in them as we get a finite 
result at zero order in Ec/Ej. In the opposite limit 
Ej <^ Ec, the situation is quite different. If one ne- 
glects Ej, the current vanishes because n is a constant 
of motion. The first non vanishing term must be of order 
Ej /Ec ■ Below we give the analytical expression for the 
current in the limit of strong dephasing (for = tR = tj 
and = = tc) including only the leading order in 
Ej/Ec- The strong dephasing Hmit refers to the con- 
dition ^jtj,Tjtj,^ctc 3> 1, which allows a series ex- 
pansion of the Josephson current at first order in e~^-'^' , 
e"^'*'^ and e-'^^^tc. 
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^strong 



- tanh I — ^ 



El 

Ec 



-ictc 



sin(2x) (cos(26')e-'^' 



cos(2x)e" 



(A4) 



r 



In this cases the Josephson energy does not any more 
enter the current through the combination Ejtj as in 
the previous case, but rather it appears through the 
ratio Ej/Ec- There is an overall suppression factor 
cx {Ej/EcY. Note that the dependence on the dynami- 
cal phase X is not the same of that in Eq. (p9)) . The ap- 



proximation Ej <^ Ec also induces a hierarchy in the de- 
phasing rates, 7}°^" - > 75°^- ~ 2{Ej/Ecf -fj. Within 
such approximation, we can neglect terms of order e~^-'*-' 
with respect to 1 (or equivalently e~'^ with respect 

to 6-^"°" * ') in the current expression, leading to 



J 



Istr 



ong 



' 7^ tanh i — ^ 

T \n 



Ej_ 

Ec 



r 



Ijtj+lctc 



sin(2x) sine; 



(A5) 



The exponential suppression due to dephasing in the map in Eq. ((52|) . The starting point of the calculation, 
Josephson contact times is reduced by the factor following the idea of Ref . Is^. is to rewrite the Eq . ([5^ . 
sm{2fi) ^ [Ej/Ecf. 



Appendix B: PHASE DEPENDENT 
CORRECTIONS TO THE CLASSICAL 
DIFFUSION 



In this Appendix we derive Eq. ([55]) for the charge 
diffusion coefficient obtained for the modified Chirikov 



Pt = pt-i- Ksm{et-i-[{t + l) raoA2]<j)} 
Ot = Ot-i+pu 

(Bl) 

as a time evolution equation for the probability distribu- 
tion function in the phase space, P{9,n;t) in which we 
add a further random step. It reads 



J 



P{e,n;t)= / G{e-e',n)P{9\n + Ksm{9' -[(t+l)mod2]<l>};t~l) dx, 
Jo 

where t counts the number of kicks. The kernel G defines the random step between two kicks. 



(B2) 



G{S0,n) 




{Se-n + 2TTmy 
2^ 



00 

27T ^ 



^ — (Tm'^/2^im{S6 — n) 



(B3) 



It correspond to add a diffusive term {a/2 the diffusion 
coefficient) in the differential equation for the time evo- 
lution between two kicks. In this way we have replaced 
the deterministic dynamics of Eq. (|Bip with a stochastic 
one which reproduce the effect of chaotic dynamics. At 
the end of the calculations we can take the limit a ^ 0. 

By introducing the Fourier coefficients for the proba- 
bility distribution function P{0, n; t) defined through 



PiO,n;t) 



1 



00 



(27r)2 ^ 



dpa(,*)(p)e*('"^+P"^ (B4) 



we can write the diffusion coefficient as 



D^- (n') = - hm ( 

2t^ ' 2t p^o+ V dp 



(B5) 



We have implicitly assumed that po = Q. We now rewrite 
the time evolution in Eq. (jB2[) in the Fourier space and 
then we expand in powers of \/\/K. The ffist step is 
performed by simply inserting Eq. (|B3I IB4p into Eq. (|B2p . 
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After some algebra, we find 



J2 Ji{K\p'\)e"^'^"-^P"> e- 



,(2t+l) 



./ — — oo 



(2t-l) 



(B6) 



+ 00 



./— — 00 

p' = p + m m! = m — I sgn{p') , (B8) 



where Ji{x) is the Zth Bessel function. 

If we represent the variables p and m respectively in 
the a;-axis and y-axis in the plane, Eq. (|B8[) defines a 
path in such a plane. The calculations of the diffusion 
constant through Eq. (jBSP is therefore reduced to the 



calculation of amip) along the path defined by Eq. (|B8p . 
Indeed many path can contribute and one has to sum over 
them. As it is evident from Eq. (|B5[) only terms of a™ (p) 
linear in time are important, therefore we can consider an 
even number of periods. We assume the initial condition 
am~°\p) = Sm,oSpfi+- Then, from Eq. (|B8|) the first step 
of the path is either (0,0) ^ (-1,1), or (0,0) -> (1,-1), 
or (0, 0) (0, 0). From Eq. (jB5|) we realize that also the 
final point of the path have to be at m = 0. The number 
of steps of the path correspond to the number of kicks. 
A trivial path consists in remaining at the origin. 



(0+,0)^(0+,0) 



(0+,0). 



(B9) 



2t+l 



The contribution of this path is 



{p^O+) = [Jo{K\p\)fa^^\p) 



l-2t 



Kp 



(BIO) 



which, combined with Eq. (|B5|) . gives D ~ K^/A. This is 
the standard result corresponding to the zero order term 
in \j\fK expansion. Any other term corresponds to a 
path with steps moving away from the origin. Any such 
step in the maps in Eq. (jB6[ IB7|) involves a Bessel func- 
tion J{x ~ K), which, in the limit K \ decays like 
~ \/^/K. Therefore a perturbation in \/K is a pertur- 
bation in the number of steps away from the origin in the 
path. The path in Eq. (jB9p is the only contribution at 
zero order and it is independent on the phase difference 

The first correction to the result in Eq. (jBlOP involves 
the path (0, 0) (1, -1) (0, 1) (0, 0) and its sym- 
metric with respect to the origin of axes. By identifying 



the path with its numerical value, we write the correc- 
tions to the diffusion constant as 5D = Vi +V2, with 




(Bll) 

In calculating the contribution of the first of the two 
paths we have to consider that the steps out of the origin 
can start at any of the {2t — 2) intermediate steps 



(p ^ 0+) = (2t - 2) [Jo{Kp)f'' J^,{Kp)e-'t'M{p + l)K)e--I^J^,{Kp)e-'t'e-^l^a^^\p) 
« (2t-2) e-2»^-V2(/v). 



(B12) 



The same procedure leads to V2 = 'Pi(</' ^ ~<t>)- We can 
now determine the correction to the diffusion coefficient: 

I? = ^[l-2cos(20)J2(i^)]+O(;|), (B13) 

If we observe that we measured time in number of kicks, 
i.e. in units of T/2, Eq. (jB13p exactly coincides with 
Eq. (j53p . in which time is measured in units of T. 



Appendix C: FULL COUNTING STATISTICS OF 
COOPER PAIR SHUTTLING 

In this appendix we briefly review the technique of Full 
Counting Statistics (FCS) introduced in the seminal pa- 
per of Lcvitov and Lesovik^ and developed by Nazarov 
et al. (see Ref. l45h . We use such formalism to obtain 
Eq. (|56l[F7|) presented in the text. We present the general 
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results for charge counting statistics without considering 
spin effects^, which do not play any role in the case of 
Cooper pair transport analyzed in the paper. 

Consider a conductor in which the charge dynamics is 
determined by the Hamiltonian ifgys- We are interested 
in the statistics of the current operator / through a given 
Section of the conductor. 

We introduce the function 3(A,A, r) = 
exp[— ©(A, A, t)], we will refer to as full counting 
statistics, formally defined through 

3(A, A,t) = Tr to+A(T,0)p(0)Zil;,(T,0)l , (CI) 



p(^)(r) 



(A±A/2)/. 



(C2) 
(C3) 



It this expression the density matrix of the system is let 
evolve in time according to two different Hamiltonians for 
bra and ket. The trace of this "modified" density matrix 
is the full counting statistics. The difference between 
the two Hamiltonians is due to the presence of A ^ 0, 
that is called counting field. It can be seen by direct 
verification that 3(A, A, r) is the generating function of 
current moments: 



5r3(A,A,T) 



A,A=0 



dh-dt„ (i{ti) 



■i{tn) 

(C4) 



Let us also observe that the derivative of S(A, A, r) gives 
us the cumulants of the current instead of its moments: 



ar6(A,A,r)U^A^„ = -9;:(ln3(A,A,T))U^A^„ 



dti...dtn(i{ti)...i{tn) 



where 



(C5) 



indicates the cumulant of the distribution 



function. The FCS can be adapted to determine vari- 
ous integrated correlation functions differing in the time 
ordering of the current operators and to obtain the cor- 
relation functions at arbitrary times rather than the in- 
tegrated ones, thus revealing as a powerful tools in the 
investigation of the transport properties of a system. 

Since the FCS can provide us information about the 
transport properties of the system, it would be mean- 
ingful to have an interpretation of 3(A, A, r) directly in 
terms of charge transfer. As we are interested in the 
transport properties of a given system contacted to two 
electron reservoirs,— the properties of the system are 
fully characterized by the knowledge of the probability 
Pt{N) of transferring N charges in a given time inter- 
val T. Pt{N) is equivalcntly defined by its generating 
function x(A) = (A^) exp(iiVA). In particular all 

moments of the current distribution function can be ob- 
tained as 



dh... dt^ (/(ti) . . . = Y^ieNTPriN) ^ (-ze)™5r xWW^^ 

"* N 



(C6) 



r 



Due to the relation between the current moments and 
the function x(A), it would appear natural to interpret 
the 3(A, A, r) just constructed as the generating function 
of the probabilities Pr{N). It can be shown that 



Pr(A,iV) 



1 

2^ 



dA3(A,A,r)e^'«^ (C7) 



is a positive defined probability if 3(A, A,t), is indepen- 
dent on A, in which case also -Pt(A, TV) is. 

Despite the previous analysis there are cases in which 
an interpretation of the transport properties in terms of 
classical probabilities is not possible. Coherent charge 
transfer between superconductors is a paradigmatic ex- 
ample. The interpretation of the FCS in this case has 
been discussed in literature (see Ref. 113 and references 
therein). Here we just notice that the appearance of non- 
positive or imaginary values of Pt (A, TV) are signature of 
the coherent nature of charge transfer between supercon- 



ductors. 

We need to adapt the procedure to the case of a tunnel 
junction where the current operator is defined in terms 
of creation and annihilation operator in the two different 
leads. Due to the linear coupling between the current and 
the counting field A in Eq. (jC3p . charge conservation en- 
sures us that the counting field fulfills all properties of a 
U{1) gauge field. In other words it means that the count- 
ing field plays the same role of the electromagnetic field, 
A. In the case of a tunnel junction between two elec- 
trons reservoirs we arc interested in the current through 
the region between the two leads. The Hamiltonian is 

H= Ekc\k^A,k+HT + H^T^ (C8) 



A=L,R 



Ht = Y ThJiC^R,kCL.h , 



(C9) 



k,h 



where Ek is the fc-dependcnt energy of electrons in the 
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reservoirs, T^.^, the tunneling amplitudes and CL^n^^k, 
'^L{Fi) k respectively the annihilation and creation op- 
erator for electrons in the left (right) reservoir. As the 
Hamiltonian is written in terms of operator defined at 
spatially separated points -that is at left (L) at right 
leads (i?)-, we have to insert the counting field in the 
only way compatible with the U{1) symmetry of the the- 
ory. It means that the counting field must enter through 
the Wilson line^, 

exp (^ij A dr^ ^ e-'(^+^/2) . (CIO) 

It results in a A-dependent Hamiltonian i/A+A/2 obtained 
by the Hamiltonian in Eq. (jCSp after the replacement 

Ht ^ e-*(^+^/2)^^ (Cll) 

Indeed the integral in Eq. (jClOp depends on the specific 
path going from the left to the right lead, but, in our case, 
the parameter A ± A/2 is defined in terms of the whole 
integral and, thus, the path-dependence in Eq. (jClOp is 
irrelevant. We can now construct the FCS for the Cooper 
pair shuttle. 

The model system we refer to is defined in Section llll 
For sake of simplicity, here we consider only the case 
tL=tR = tj, Ef^ = Ef^ = Ej, and = = tc- 
We consider the counting statistics for electrons passing 
through a counter at left lead. The charge transfer in the 



Cooper pair shuttle takes place through the Josephson 
effect, i.e. coherent tunneling of Cooper pairs. Accord- 
ing to the general procedure just described, we have to 
construct i?A.+ A/2 fo r tunneling of Cooper pairs by 
modifying Eq. (|C11[) . In the Hamiltonian of the Cooper 
pair shuttle, Eq. ([47|) . the tunneling term corresponds 
to the creation of a Cooper pair (of charge 2e) into the 
grain, therefore the modification of the Hamiltonian in 
Eq (jClip is given by 

\n + 1) {n\ ^ e-''(2A+A) |„ + ^ (ci2) 

or equivalcntly 

i?A±^/2 = M^L ^ 0L + 2A ± A) . (C13) 

In this way A is counting the charge in units of e. Note 
that A can be reabsorbed into </)l by the redefinition 
<j>L + A ^ (f>L. Therefore, as long as the FCS depends on 
0, it will depend on A and the interpretation in terms of 
probabilities that a certain number of Cooper pairs have 
traversed the shuttle will be impossible. The construc- 
tion of the FCS including counting fields both at left and 
right contacts is straightforward and gives 

H\L/2.Xn/2 = Ho{(j)L 4>L - >^L, 4>B. ^ 4>B^ ~ A_r) , 

(C14) 

as presented in Eq. ([57)) . 
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If the Cooper pair shuttling is achieved by means of time- 
dependent fiuxes, as described in Subsection III Al Ec(t) is 
a constant during the whole cycle. 

The high temperature limit is intended Tj, ^ Ej, but, in 
any case, Tj <^ Ec, such that we can limit our discussion 
to two charge states. 

A Gaussian distribution function for the variables Ati,(j) 
cannot be assumed; it is inconsistent with the previous 
condition, 'P{Atb{i)) = for Atb{i) < 0, and, in fact, gives 
rise to incurable divergences due to the presence of the 
term 1/T in the expressions to be averaged. 
Periodically driven Josephson junction have been already 
suggested to study quantum chaos in Ref. |49| . 
"'"^ The generalization to the multi-terminal case is straight- 
forward. 



